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Multifractal  characteristics  of  the  upper  layer  (above  140  m  depth)  thermal  structure  in  the  western  Philippine  Sea  near 
Taiwan  are  analyzed  using  high-resolution,  digital  thermistor  chain  data.  The  power  spectra  at  all  the  depths  have  multi¬ 
scale  characteristics  with  the  spectral  exponent  p  in  the  range  of  (1,  2),  which  indicates  nonstationary  with  stationary 
increments.  The  graph  dimension  varies  from  higher  values  such  as  1.71  (in  sublayer:  60  m),  to  lower  values  such  as  1.59 
(in  second  themiocline:  120  m).  However,  the  information  dimension  varies  slightly  from  0.929  to  0.941.  The  multi¬ 
dimensional  structure  is  stronger  in  nonstationarity  (graph  dimension)  than  in  intermittency  (information  dimension).  These 
results  provide  useful  information  for  investigating  turbulence  structure  in  the  upper  layer  of  the  western  Philippine  Sea. 

[Key  words:  Multifractal  structure,  high-order  structure  function,  singular  measure,  power  law,  multifractal  plane] 


1.  Introduction 

The  western  Philippine  Sea  (WPS)  near  eastern 
Taiwan  (Fig.  1)  is  the  region  that  the  Philippines 
Current  flows  northward  from  just  north  of  Mindanao 
(the  southernmost  of  the  main  island  groups  of  the 
Philippines)  at  around  (10°N,  128°E)  to  Taiwan  where 
it  continues  as  the  Kuroshio  -  the  most  evident 
western  boundary  current  in  the  North  Pacific  Ocean. 
After  becoming  the  Kuroshio,  it  makes  a  slight 
excursion  into  the  South  China  Sea  through  the  Luzon 
Strait,  continues  northward,  and  is  fed  from  the  North 
Equatorial  Current,  from  which  it  continues  to  entrain 
water  along  its  way.  The  Kuroshio  current  system  has 
various  stages  of  modification  that  transports  warm 
water  to  the  North  Pacific.  Because  of  the  import  and 
modification  of  water  masses  a  large  number  of 
regional  water  types  can  be  encountered. 

Eddy  activity  is  often  observed  in  WPS  near 
Taiwan.  For  example,  a  large  warm  eddy  lies  east  of 
the  Luzon  Strait,  a  seemingly  permanent  feature1 3 . 
Although  the  WPS  is  the  only  important  source  of 
open-ocean  water  to  the  South  China  Sea,  the 
hydrographic  and  chemical  properties  of  the  waters  on 
both  side  of  the  Luzon  Strait  are  quite  different1. 
Recent  observations  also  show  no  evidence  of 
branching  of  Kuroshio  into  the  South  China  Sea, 
although  Kuroshio  water  is  often  found  inside  the 
South  China  Sea4.  It  is  likely  that  intrusion  of  the 
WPS  water  in  the  upper  600  m  is  caused  by  the 
(density)  frontal  instability  of  Kuroshio5. 


Oceanic  density  fronts  often  form  cusp-like 
corrugation  at  the  sea  surface.  These  frontal  cusps, 
ranging  in  scale  from  meters  to  kilometers,  are  often  a 
prominent  feature  in  radar  images  of  coastal  fronts. 
The  complex  three-dimensional  flow  structure 
underlying  the  surface  cusp  has  been  difficult  to  study 
despite  many  past  observations.  Two  major  features, 
nonstationarity  and  intermittency,  should  be  first 
investigated  in  order  to  understand  oceanic  frontal 
instability.  Question  arises:  How  can  we  determine 


Fig.  1 — Topography  of  the  western  Philippine  Sea  and 
surrounding  areas.  The  inserted  box  shows  the  track  of  CMB. 
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upper  ocean  nonstationarity  and  intermittency  from 
temperature  observations?  This  paper  uses  a 
multifractal  analysis  on  a  high-resolution  temperature 
dataset  to  obtain  the  nonstationarity  and  intermittency 
of  the  upper  layer  (140  m  deep)  in  the  WPS. 

2.  Coastal  Monitoring  Buoy  with  Attached 
Thermistors 

The  coastal  monitoring  buoy  (CMB)  was  deployed 
by  the  U.S.  Naval  Oceanographic  Office  in  WPS 
(Fig.  2)  during  28  July  -  7  August  2005.  The  original 
objective  was  to  collect  the  data  every  10  min  near  the 
air-sea  interface.  Above  the  ocean  surface,  the  surface 
winds,  air  temperature,  and  air  pressure  were 
measured.  Below  the  ocean  surface,  the  temperature 
was  observed  at  1,3,  5,  10,  15,  and  20  m.  Fifteen 
thermistors  were  attached  to  a  wire  rope  extending 
from  the  code  of  CMB  (20  m  deep)  to  140  m  with 
high  frequency  sampling  rate  (every  15  sec)  (Fig.  3). 
The  purpose  was  to  observe  detail  of  the  upper  layer 
structure.  Figure  4  shows  contour  plot  of  temperature 
ranging  from  2 1  °C  to  29°C  with  1 ,0°C  increment  on  a 
vertical  cross-section  along  the  track  (Fig.  2).  The 
total  length  of  the  cross-section  is  229.14  km.  The 
CMB  travels  3.82  m  every  15  seconds. 

The  temperature  shows  a  multi-scale  variability 
with  high  irregularity.  The  surface  mixed  layer  is  very 


Costal  Monitoring  Buoy  Track 


Fig.  2 — Track  of  CMB  (from  July  28  to  August  7,  2005)  deployed 
by  the  Naval  Oceanographic  Office. 


thin  (depth  around  5  m)  on  28  July  2005  with 
temperature  about  28.5  °C.  Below  the  surface  mixed 
layer,  two  thermoclines  appear  with  the  first 
thermocline  to  the  depth  of  50  m.  The  vertical 
gradient  in  the  first  thermocline  is  around  0.1°C/m.  A 
relatively  uniform  sublayer  (24  °C)  exists  below  the 
first  thermocline  from  50  to  130  m.  Below  the 
uniform  sublayer,  there  is  a  second  thermocline  with  a 
vertical  temperature  gradient  around  0.04  nC/m.  As 
time  approaches,  the  surface  mixed  layer  deepens  and 
the  first  thermocline  descends.  These  processes 


Fig.  3 — The  coastal  monitoring  buoy  used  in  the  WPS  survey. 
Fifteen  thermistors  are  attached  to  a  wire  rope  extending  from  the 
code  of  CMB  (20  m  deep)  to  140  m  with  high  frequency  sampling 
rate  (every  15  sec). 
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contain  small  scale  fluctuations.  The  surface  mixed 
layer  reaches  70  m  on  5  August  2005,  when  an 
evident  cooling  occurs  with  the  mixed  layer 
temperature  reducing  to  25.5  °C.  At  the  same  time,  the 
two  thermoclines  merge  to  a  single  themocline  with 
the  vertical  gradient  of  0.06  °C/m. 

3.  Power  Spectra 

Determination  of  the  stationarity  is  the  first  step 
towards  understanding  the  inherent  thermal  variability 
and  statistical  properties  identified  from  this  high- 
resolution  temperature  data  with  multi-layer 
structures.  For  a  given  depth,  the  temperature  data  is  a 
function  of  the  horizontal  coordinates,  x, 

Tt  =T(xi),  x,=il,  i  =  0,1,..., A,  A  =L/l,  ...  (1) 

where  /  =  3.8  m,  is  the  horizontal  resolution  of  the 
data,  and  L  is  the  total  horizontal  scale  of  the  data  set. 
Our  data  set  has  60,300  temperature  profiles  (A  = 
60,300). 


Spectral  analyses  of  temperature  field, 

EJ  =  E(kj ),  kj  =  j/ L,  j  =  1,2,..., A/2  ...  (2) 

at  all  depths  were  conducted.  Since  large  number  of 
input  values  (60,  300)  is  considered,  the  Bartlett 
window  is  used  after  the  Fast  Fourier  Transform  in 
order  to  reduce  the  spectral  leakage  in  the 
wavenumber  domain.  For  the  sake  of  brevity  and  to 
elucidate  the  important  points,  only  spectra  at  six 
depths  are  shown  in  Fig.  5. 

For  a  scaling  process,  one  expects  power  law 
behavior6,7, 

£(*)«*-',  ...(3) 

over  a  large  range  of  wavenumber  k.  The  spectral 
exponent  p  contains  information  about  the  degree  of 
stationarity  of  the  data.  If  P  <  1 ,  the  field  is  stationary; 
if  1  <  P  <  3,  the  field  contains  nonstationary  signal 
with  stationary  increments  and  in  particular,  the 


Temperature  Contoured  all  °C 
Beginning  28-Jul-2005  09  08  00 


07/28  07/23  07/30  07/31  08/01  00/02  08/03  08/04  08/05  08/06  08/07 


Fig.  4 — Temperature  cross  section  obtained  from  the  CMB  data  collected  along  the  track  shown  in  Fig.  2. 


lofl2  E(k)  log2  E(k)  loflj  E(k) 
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small-scale  gradient  field  is  stationary;  if  j3  >  3,  the 
field  is  nonstationary  with  nonstationary  increments. 
The  power  spectra  for  all  the  depths  (Fig.  5)  have 
multi-scale  characteristics  with  the  spectral  exponent 
j8  in  the  range  of 

l<(3<2,  ...  (4) 

which  means  the  temperature  field  of  the  WPS  upper 
layer  is  nonstationary  with  stationary  increments. 

4.  Structure  Functions 

Since  the  thermistor  chain  data  set  has  stationary 
increments,  we  should  study  the  statistical 
characteristics  of  the  gradient  field, 

I  AF(r ; x)  1=1  T (xi+r )-T(xj)\,  i  =  0,1,...,A -r  ,  ...  (5) 

where  r  denotes  the  lag  between  two  data  points.  The 
gradient  may  bound  to  create  certain  amount  of  noise 
and  with  different  lags.  However,  the  noise  constitutes 
the  high-frequency  portion  of  the  power  law  relation, 
it  usually  does  not  change  the  structure  of  the  power 
law.  Obviously,  r  is  inversely  proportional  to  the 
wavenumber  k , 

1 

r  oe  —  . 
k 

Structure  functions  are  used  to  depict  the  thermal 
variability.  The  ^th-order  structure  function  is  defined 
by  the  mean  of  the  qth  power  of  the  gradient 
field  I  AT(r,  x)  I, 

S(r,<?)^{|Air(r;x)P)  =  ~1—  ^|AT(r;x)f.  ...  (6) 

'  '  A  —  r  i=0 

For  example,  r  =  1 .  q  =  1 ,  the  structure  function 
5(1, 1 )  =  — X |r (xM ) - 7X*i )| ,  ...  (7) 

A  —  1  j=o 

represents  the  average  magnitude  of  gradient. 

Near-linear  dependence  of  Iog2[S(r,^)]  on  log2(r) 
is  found  with  different  ^-values  from  0.5  to  4.0 
(Fig.  6).  The  straight  lines  with  different  slopes  show 
that  the  structure  functions  for  the  upper  ocean 
temperature  in  the  WPS  satisfies  the  power  law 


S(r,q)ocriw,  ...  (8) 

with  the  exponent  £(q)  depending  on  q.  Since  S(r,  0) 
=  1,  the  power  £(q)  should  be  0  for  q  =  0.  Our 
computation  (Fig.  6)  agrees  quite  well  with  earlier 
studies  in  the  southwestern  Greenland  Sea,  Iceland 
Sea,  and  Norwegian  Sea6:  £  (q)  is  monotonically  and 
near-linearly  increasing  with  q.  Thus,  it  may  be 
represented  by 

£(q)  =  H(q)q,  ...  (9) 

where  H(q)  is  nearly  a  constant. 

The  power  of  the  structure  function,  £  (<7) ,  is 
computed  for  the  whole  dataset.  Figure  7  shows  the 
dependence  of  £  (q)  on  q  for  selected  levels.  All  the 
curves  (near-linear)  converges  at  q  =  0  [£(0)  =  0], 
and  show  increase  of  C,  (q)  with  q  for  all  depths. 
Such  an  increase  is  quite  complicated  as  q  <  3  and 
shows  three  different  patterns  as  q  >  3:  (a)  slow 
increasing  with  q  in  the  sub-layer  (50  -  80  m),  (b) 
intermediate-rate  increasing  with  q  in  the  upper 
thermocline  (35  -  45  m)  and  lower  thermocline  (100  - 
140  m),  (c)  fast  increasing  with  q  in  the  mixed  layer 
(25 -30  m). 

The  structure  function  for  q  ~  1,  S(r,  1),  is  often 
used  to  determine  the  statistical  characteristics  of  the 
data  such  as  stochastically  continuous  and 
stationarity.  When  the  structure  function  S(r,  1 )  does 
not  depend  on  r 

S(r,l)  =  const,  ...  (10) 

which  denotes  exact  stationary.  Eq.  (10)  is  equivalent 

to 

f(l)  =  //(l)  =  0. 

Thus,  M\  =  //( 1 ),  represents  statistical  characteristics 
of  data.  Since  the  power  £(1)  [i.e.,  //(!)]  varies  with 
depth  (Fig.  7)  between  0.29  at  60  m  depth  (sub-layer) 
and  0.41  at  140  m  depth  (second  thermocline),  the 
thermal  field  in  the  WPS  is  nearly  stationary. 

If  g(T)  is  defined  as  the  graph  of  T(x),  whose 
dimension  can  be  defined8  by 

®s(T)  ~2  —  Hl.  ...(11) 


INDIAN  J.  MAR.  SCI.,  VOL.  36,  NO.  2,  JUNE  2007 


4  5  6  7  8  9  10  11  12  13  4  5  6  7  8  9  10  11  12  13 


tog2(r)  k>g2(r) 


Fig.  6 — Structure  functions  for  different  (/-values  at  depths  of  (A)  25  m,  (B)  40  m,  (C)  50  m,  (D)  75  m,  (E)  100  m,  and  (F)  140  m. 
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Dependence  of  the  Structure  Function  Zeta(q.nlev) 


q 


Fig.  7 — Dependence  of  the  structure  function’s  power.  C  (q) ,  on  q  and  depth. 


For  a  stochastically  continuous  (H\  =  1)  data  T(x),  the 
graph  is  a  smooth  curve,  whose  dimension,  Dg(T)  =  1. 
If  the  graphs  g(T)  fill  the  whole  space  (exact 
stationary), 

Dgm  =  2,  ...  (12) 

which  corresponds  to  H \  =  0.  The  graph  dimension  of 
upper  thermal  field  of  the  WPS  varies  (multi¬ 
dimension  structure)  from  higher  values  such  as  1 .7 1 
(sublayer:  60  m),  to  lower  values  such  as  1.59  (second 
thermocline:  120  m). 

5.  Singular  Measures 

Normalized  small  scale  absolute  gradient 


is  used  to  identify  the  intermittency  of  the  thermal 
field.  The  running  average  of  r  normalized  values  are 
computed  by 

i  i+r-1 

e(r-x.)  =  -  y  £(1;*.),  i  =  0,1,...,  A  -  r  .  ...  (14) 

r  y=. 

The  mean  of  the  q\h  power  of  e(r;jci  ) 

M(r,?)  =  (£(r;r,),)  =  -!-^[£(r;i:i)J'  ,  ...  (15) 

'  A -r£o 

is  defined  as  the  c/th-order  singular  measure. 
Obviously,  for  q  =  0, 


«(!;*,) 


lA7-(l;x,.)| 

(|A7’(1;JC[.)|)’ 


{|AT(l;^)|)  =  |X|A7’(lu1.)|, 

A  i= o 


M(r,  0)  =  1 .  ...(16) 

For  q  =  1 , 

M  (r,l)  =  (e(r;x:)}  =  —^—  ^[fi(ru,)] 

A  —  r  " 


...  (13) 


loflj  (M(r,q)l  log2  (M(r,q)]  tog2{M(r,q)I 
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Fig.  8 — Singular  measures  for  different  ^-values  at  depths  of  (A)  25  m,  Bb)  40  m,  (C)  50  m,  (D)  75  m,  (E)  100  m,  and  (F)  140  m. 
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A -r3 


i  i+r-1 


=  1. 


...  (17) 


d2K{q) 
dq 2 


>0, 


...  (20) 


The  singular  measures  are  computed  for  all  depths. 
For  simplicity,  M(r,  q)  for  six  depths  (25,  40,  50,  75, 
100,  140  m)  is  given  here  (Fig.  8).  Near-linear 
dependence  of  iog2[M(r,q)]  on  log2(r)  is  found 

with  different  ^-values  from  0.5  to  4.0.  The  straight 
lines  with  different  slopes  show  that  the  singular 
measures  with  various  q  for  the  upper  layer 
temperature  in  the  WPS  satisfies  the  power  law 

M(r,q)  r~K{9\  q>  0,  ...  (18) 


for  all  q  and 

K{q)<  0  only  if  0  <  q  <  1  ...  (21) 

which  reflects  the  fact  that,  in  this  range,  taking  a  qth 
power  necessarily  reduces  the  fluctuation  of  £  (r;  jc,); 
and  otherwise 

K(q)>  0,  if  q>0.  ...  (22) 

Following  Chu6,  we  may  define  a  function 


C{q)  =  ^-.  ...  (23) 

q- 1 

For  q  — » 1,  we  use  I'Hosptal'a  rule  to  define  a 
straightforward  measure  of  inhomogeneity  in  the 
sense  of  singular  measure 

0.6 

i 

0.5  1 

0,4  ( 

5*0.3 

sr 

0.2 

0.1 

0 

0  0.5  1  1.5  2  2.5  3  3.5  4 

q 


with  the  power  K(q)  varying  with  q.  From  Eqs.  (16) 
and  (17),  we  have 


AT(0)  =  ^(1)  =  0. 


(19) 


Several  characteristics  are  found  from  Fig.  9:  The 
power  K(q)  is  a  convex  function 


Fig.  9 — Dependence  of  the  singular  measure’s  power,  K{q) ,  on  q  and  depth. 
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Mean  multi-fractal  Plane  (H1,C1) 


Fig.  10 — Mean  multi-fractal  plane  (Hu  CO  for  the  western  Philippine  Sea  thermal  field. 


C,  =C(1)  =  /T(1)^0,  ...(24) 

which  is  called  the  intermittency  parameter.  The 
larger  the  value  of  C\,  the  lager  the  intermittency  and 
singularity  the  data  set  has.  The  intermittency 
parameter  Ci  varies  from  0.059  to  0.07 1 . 

6.  Mean  Multifractal  Plane 

Statistical  moments  of  all  orders  (t/th  power)  can  be 
computed,  however,  the  first  moment  (q  =  1)  provides 
important  information  about  the  data.  Roughness 
indicates  differentiability  while  singularity  indicates 
uniqueness.  Sparseness  may  be  inferred  as  either  lack 
of  information  or  no  continuity,  but  it  does  not  mean 
singularity.  Roughness  constitutes  with  information 
and  shows  continuity  and  for  that  reason  it  can  be 
differentiable.  The  parameter  C j  measures  the  degree 
of  intermittency  in  the  system,  while  Hi  measures  its 
degree  of  non  stationary. 

The  plot  of  C\  versus  Hu  called  the  mean 
multifractal  plane,  shows  the  degree  of  nonstaionarity 
and  intermittency  (Fig.  10).  Both  parameters  have 
geometrical  interpretations  as  co-dimensions: 


information  dimension  and  graph  dimension6.  The 
information  dimension,  represented  by  (1  -  C\),  is  a 
first-order  estimate  of  sparseness  of  strong  gradient 
distributed  in  the  system.  The  graph  dimension, 
represented  by  (2  -  H\),  is  a  first-order  estimate  of 
roughness  in  the  system.  Both  parameters  have 
analytical  meanings:  C\  (sparseness)  is  related  to 
singularity  and  H\  (roughness)  is  connected  to  the 
lack  of  differentiability.  The  WPS  upper  layer  thermal 
field  shows  multifractal  characteristics  (Fig.  10)  in 
terms  of  C\  and  H\.  We  find  larger  variation  in  //, 
(0.29-0.41)  than  in  C,  (0.059-0.071).  Identification  of 
(Hu  Ci)  ranges  helps  to  select  realistic  model  to 
describe  the  field. 

7.  Conclusion 

The  upper  layer  thermal  structure  in  the  western 
Philippine  Sea  has  the  following  multifractal  features: 

(1)  The  power  spectra  at  all  the  depths  have  multi¬ 
scale  characteristics  with  the  spectral  exponent  /?  in 
the  range  of  (1,  2).  This  means  the  temperature  field  is 
nonstationary  with  stationary  increments. 
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(2)  The  structure  function  has  multifractal 
characteristics.  The  powers  of  the  4th -order  structure 
function  show  three  different  patterns:  (a)  slow 
increasing  with  q  in  the  sub-layer  (50  -  80  m), 
(b)  intermediate-rate  increasing  with  q  in  the  upper 
thermocline  (35  -  45  m)  and  lower  thermocline  ( 100  - 
140  m),  (c)  fast  increasing  with  q  in  the  mixed  layer 
(25  -  30  m). 

(3)  Two  parameters  for  the  first  moment  (q  =  1) 
structure  function  and  sigular  measure,  H\  and  C are 
important  to  represent  the  multifractal  characteristics 
in  terms  of  co-dimensions  (information  and  graph). 
The  graph  dimension,  represented  by  (2  -  H ,),  is  a 
first-order  estimate  of  roughness  in  the  system.  The 
information  dimension,  represented  by  (1  -  C|),  is  a 
first-order  estimate  of  sparseness  of  strong  gradient 
distributed  in  the  system.  Both  parameters  have 
analytical  meanings:  C)  (sparseness)  is  related  to 
singularity  and  H 1  (roughness)  is  connected  to  the 
lack  of  differentiability. 

The  graph  dimension  of  the  upper  thermal  field  of 
the  western  Philippine  Sea  varies  (multi-dimension 
structure)  from  higher  values  such  as  1.71  (sublayer: 
60  m),  to  lower  values  such  as  1.59  (second 
thermocline:  120  m).  However,  the  information 
dimension  of  the  western  Philippine  Sea  upper 
thermal  field  varies  slightly  from  0.929  to  0.941.  This 
indicates  that  singularity  is  very  low. 
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